home *** CD-ROM | disk | FTP | other *** search
/ The Very Best of Atari Inside / The Very Best of Atari Inside 1.iso / mint / mntlb20 / lib / _muldf3.cpp < prev    next >
Text File  |  1992-03-30  |  4KB  |  194 lines
  1. | mjr: not needed on the TT
  2.  
  3. #ifndef    __M68881__
  4. # ifdef    sfp004
  5.  
  6. | double precision floating point stuff for Atari-gcc using the SFP004
  7. | developed with gas
  8. |
  9. | double precision multiplication
  10. |
  11. | M. Ritzert (mjr at dfg.dbp.de)
  12. |
  13. | 4.10.1990
  14. |
  15. | no NAN checking implemented since the 68881 treats this situation "correctly",
  16. | i.e. according to IEEE
  17.  
  18. | addresses of the 68881 data port. This choice is fastest when much data is
  19. | transferred between the two processors.
  20.  
  21. comm =     -6
  22. resp =    -16
  23. zahl =      0
  24.  
  25. | waiting loop ...
  26. |
  27. | wait:
  28. | ww:    cmpiw    #0x8900,a0@(resp)
  29. |     beq    ww
  30. | is coded directly by
  31. |    .long    0x0c688900, 0xfff067f8
  32.  
  33.     .text
  34.     .even
  35.     .globl    __muldf3, ___muldf3
  36.  
  37. __muldf3:
  38. ___muldf3:
  39.     lea    0xfffa50,a0
  40.     movew    #0x5400,a0@(comm)    | load first argument to fp0
  41.     cmpiw    #0x8900,a0@(resp)    | check
  42.     movel    a7@(4),a0@
  43.     movel    a7@(8),a0@
  44.     movew    #0x5423,a0@(comm)
  45.     .long    0x0c688900, 0xfff067f8
  46.     movel    a7@(12),a0@
  47.     movel    a7@(16),a0@
  48.     movew    #0x7400,a0@(comm)    | result to d0/d1
  49.     .long    0x0c688900, 0xfff067f8
  50.     movel    a0@,d0
  51.     movel    a0@,d1
  52.     rts
  53.  
  54. # else    sfp004
  55.  
  56. | double floating point multiplication routine
  57. |
  58. | written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet).
  59. | Based on a 80x86 floating point packet from comp.os.minix, written by P.Housel
  60. |
  61. |
  62. | Revision 1.2, kub 01-90 :
  63. | added support for denormalized numbers
  64. |
  65. | Revision 1.1, kub 12-89 :
  66. | Ported over to 68k assembler
  67. |
  68. | Revision 1.0:
  69. | original 8088 code from P.S.Housel
  70.  
  71. BIAS8    =    0x3FF-1
  72.  
  73.     .text
  74.     .even
  75.     .globl    __muldf3, ___muldf3
  76.  
  77. __muldf3:
  78. ___muldf3:
  79.     lea    sp@(4),a0
  80.     moveml    d2-d7,sp@-
  81.     moveml    a0@,d4-d5/d6-d7 | d4-d5 = v, d6-d7 = u
  82.     subw    #16,sp        | multiplication accumulator
  83.  
  84.     movel    d6,d0        | d0 = u.exp
  85.     swap    d0
  86.     movew    d0,d2        | d2 = u.sign
  87.     lsrw    #4,d0
  88.     andw    #0x07ff,d0    | kill sign bit
  89.  
  90.     movel    d4,d1        | d1 = v.exp
  91.     swap    d1
  92.     eorw    d1,d2        | d2 = u.sign ^ v.sign (in bit 31)
  93.     lsrw    #4,d1
  94.     andw    #0x07ff,d1    | kill sign bit
  95.  
  96.     andl    #0x0fffff,d6    | remove exponent from u.mantissa
  97.     tstw    d0        | check for zero exponent - no leading "1"
  98.     beq    0f
  99.     orl    #0x100000,d6    | restore implied leading "1"
  100.     bra    1f
  101. 0:    addw    #1,d0        | "normalize" exponent
  102. 1:    movel    d6,d3
  103.     orl    d7,d3
  104.     beq    retz        | multiplying by zero
  105.  
  106.     andl    #0x0fffff,d4    | remove exponent from v.mantissa
  107.     tstw    d1        | check for zero exponent - no leading "1"
  108.     beq    0f
  109.     orl    #0x100000,d4    | restore implied leading "1"
  110.     bra    1f
  111. 0:    addw    #1,d1        | "normalize" exponent
  112. 1:    movel    d4,d3
  113.     orl    d5,d3
  114.     beq    retz        | multiplying by zero
  115.  
  116.     addw    d1,d0        | add exponents,
  117.     subw    #BIAS8+16-11,d0    | remove excess bias, acnt for repositioning
  118.  
  119.     clrl    sp@        | initialize 128-bit product to zero
  120.     clrl    sp@(4)
  121.     clrl    sp@(8)
  122.     clrl    sp@(12)
  123.     lea    sp@(8),a1    | accumulator pointer
  124.  
  125. | see Knuth, Seminumerical Algorithms, section 4.3. algorithm M
  126.  
  127.     swap    d2
  128.     movew    #4-1,d2
  129. 1:
  130.     movew    d5,d3
  131.     mulu    d7,d3        | mulitply with bigit from multiplier
  132.     addl    d3,a1@(4)    | store into result
  133.     movew    d4,d3
  134.     mulu    d7,d3
  135.     movel    a1@,d1        | add to result
  136.     addxl    d3,d1
  137.     movel    d1,a1@
  138.     roxlw    a1@-        | rotate carry in
  139.  
  140.     movel    d5,d3
  141.     swap    d3
  142.     mulu    d7,d3
  143.     addl    d3,a1@(4)    | add to result
  144.     movel    d4,d3
  145.     swap    d3
  146.     mulu    d7,d3
  147.     movel    a1@,d1        | add to result
  148.     addxl    d3,d1
  149.     movel    d1,a1@
  150.  
  151.     movew    d6,d7
  152.     swap    d6
  153.     swap    d7
  154.     dbra    d2,1b
  155.  
  156.     swap    d2        | [TOP 16 BITS SHOULD BE ZERO !]
  157.  
  158.     moveml    sp@(2),d4-d7    | get the 112 valid bits
  159.     clrw    d7        | (pad to 128)
  160. 2:
  161.     cmpl    #0x0000ffff,d4    | multiply (shift) until
  162.     bhi    3f        |  1 in upper 16 result bits
  163.     cmpw    #9,d0        | give up for denormalized numbers
  164.     ble    3f
  165.     swap    d4        | (we''re getting here only when multiplying
  166.     swap    d5        |  with a denormalized number; there''s an
  167.     swap    d6        |  eventual loss of 4 bits in the rounding
  168.     swap    d7        |  byte -- what a pity 8-)
  169.     movew    d5,d4
  170.     movew    d6,d5
  171.     movew    d7,d6
  172.     clrw    d7
  173.     subw    #16,d0        | decrement exponent
  174.     bra    2b
  175. 3:
  176.     movel    d6,d1        | get rounding bits
  177.     roll    #8,d1
  178.     movel    d1,d3        | see if sticky bit should be set
  179.     orl    d7,d3        | (lower 16 bits of d7 are guaranteed to be 0)
  180.     andl    #0xffffff00,d3
  181.     beq    4f
  182.     orb    #1,d1        | set "sticky bit" if any low-order set
  183. 4:    addw    #16,sp        | remove accumulator from stack
  184.     jmp    norm_df        | (result in d4/d5)
  185.  
  186. retz:    clrl    d0        | save zero as result
  187.     clrl    d1
  188.     addw    #16,sp
  189.     moveml    sp@+,d2-d7
  190.     rts            | no normalizing neccessary
  191.  
  192. # endif    sfp004
  193. #endif    __M68881__
  194.